Fate of Dirac points in a vortex superlattice 
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We consider noninteracting fermions on the honeycomb lattice in the presence of a magnetic 
vortex superlattice. It is shown that depending on the superlattice periodicity, a gap may open at 
zero energy. We derive an expression of the gap in the small-flux limit but the main qualitative 
features are found to be valid for arbitrary fluxes. This study provides an original example of a 
metal-insulator transition induced by a strongly modulated magnetic field in graphene. At the same 
time our results directly apply to Kitaev's honeycomb model in a vortex superlattice. 
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After its experimental discovery by Geim and 
Novoselov [ ] in 2004, graphene's electronic properties 
received much attention (see Rcf . [ ] for a review) . How- 
ever, the band structure of this honeycomb lattice in the 
tight-binding approximation has been known for several 
decades, [3] and its modifications in the presence of a uni- 
form magnetic field were investigated more than 20 years 
ago by Rammal [4]. One of the most salient features 
of the zero-field spectrum is the existence of a point-like 
Fermi surface at zero energy, the celebrated Dirac points, 
giving rise to a relativistic dispersion in their neighbor- 
hood (the so-called Dirac cones). Interestingly, these dis- 
crete zero-energy states are still present when a uniform 
magnetic field, is added [4, 5]. The stability of these 
states has led several groups to analyze the influence of 
a nonuniform magnetic field and it is now commonly ac- 
cepted that a smoothly modulated magnetic field is not 
sufficient to open a gap at zero energy [6—8]. 

In this paper, we show that it is actually possible 
to open this gap by considering the opposite limit of a 
strongly modulated magnetic field. In this case, unlike 
previous studies [6-8], one cannot neglect the coupling 
between Dirac cones, which is directly responsible for 
this dramatic effect. As a consequence, the opening of 
the gap does not require the simultaneous presence of 
a scalar and a vector potential. To analyze this prob- 
lem, we consider a vortex superlattice with fluxes ±0 
as depicted in Fig. 1. Our choice is motivated by the 
commensurability of the triangular and hexagonal struc- 
tures and by the fact that this alternated pattern leads to 
the smallest possible unit cell of the superlattice. In the 
small-</> limit, we show that although the system remains 
gapless at first order, a gap proportional to </> 2 may open, 
providing a nice example of a metal-insulator transition 
induced by a magnetic field in the honeycomb lattice. We 
derive the necessary and sufficient condition to open this 
gap in terms of the superlattice periodicity, and we give 
an expression of the gap at order two in the small-</> limit. 
When the size of the superlattice unit cell increases, i.e., 
in the limit of vanishing vortex density v, we find that 



FIG. 1. (Color online) A piece of the £(p — l,q — 1) magnetic 
vortex superlattice spanned by primitive vectors foi and b-2. 
Vectors a\ and ai are primitive vectors of the bare honeycomb 
lattice. Light-center (green) and dark-center (red) plaquettes 
contain a flux +</> and — 0, respectively, whereas white pla- 
quettes are flux-free. Blue links with arrows indicate oriented 
hopping terms "carrying" the flux. 



the gap vanishes as ulnv^ 1 . 

Although obtained in a perturbative framework, our 
conclusions remain qualitatively valid for arbitrary fluxes 
as checked by exact diagonalizations. Furthermore, at 
4> = ^o/2 (0o being the elementary flux quantum), the 
same gap-opening mechanism applies to the celebrated 
Kitaev model [9] studied in the context of topologically 
ordered systems. 

The starting point of our study is the following tight- 
binding Hamiltonian 



ff=-X>ji*>oi, 



(i) 



(i,j) 



where \i) denotes a spinless-electron state localized on 
site i. The sum is performed over all nearest-neighbor 
sites of the honeycomb lattice and the hopping term in 
the presence of a vector potential A is given by the so- 
called Peier Is substitution [ I ]: tij — te^o Ji ' . Thus, 
setting the flux and energy scales to unity (cpo = t = 1), 
the (oriented) product of the hopping terms over a closed 



loop is simply e 2l7r< ^ where is the dimensionless magnetic 
flux inside the corresponding loop. 

The vortex superlattice considered here is defined as 
follows. Let us assume that there is a flux +<fi in the ele- 
mentary plaquette centered in r. Then, the superlattice 
C{p, q) is generated by requiring that the plaquette lo- 
cated at r + 6i/2 contains a flux — and the one located 
at r + &2 contains a flux +</>, where 



61 = 2(poi +qa 2 ), 

b 2 = -<?cti + {p + q)a 2 . 



(2) 
(3) 



Vectors a\ and Gi2 are primitive vectors of the honeycomb 
lattice (see Fig. 1 for the case p = q = 1) and (p, q) are 
positive integers. In the following, without loss of gener- 
ality, we only consider the case p ^ q. It is straightfor- 
ward to check that the total flux per unit cell of £(p, q), 
spanned by b\ and b 2l is zero. In addition, the vortex 
density defined as the number of vortices per unit cell is 
simply given by 



(4) 



p 2 + pq + <r 



A convenient gauge choice realizing such a flux pattern 
can be obtained starting from an initial +0 plaquette 



center and by choosing ti 



,2i7rd 



for all links crossed 



by going p times in direction a\ and then q times in 
the direction a 2 . The orientation of the first link fixes 
all others since we wish to have a flux —<f> in the final 
plaquette and zero in all intermediate ones. In other 
words, one creates a string of links carrying the flux which 
connects a vortex to an antivortex. As a side remark, let 
us note that with this gauge choice, one can study any 
value of the flux without changing the size of the unit 
cell, in deep contrast with the uniform field problem. 

As for any bipartite lattice, the spectrum of H is sym- 
metric with respect to the energy e = for all <f>. For 
<j) = 0, it consists of two symmetric bands [ 



e±(k) — ±< 3 + 2cos(fc.ai) + 2cos(fc.a,2) 

I 1/2 
+2cos[fe.(ai-a 2 )]j . (5) 



These symmetric bands touch at e = when k coin- 
cides with the so-called Dirac points K = ^a\ + |a^ and 
K' = |a* + \cl 2i where a* and a 2 are primitive vec- 
tors of the reciprocal lattice associated to cti and a 2 
{a*.a,j = 2n5ij). Consequently, the energy e — is four- 
fold degenerate for = 0. Our goal is to determine the 
fate of these zero-energy states for <f> 7^ 0. 

To address this problem, we shall analyze perturba- 
tively the small-</> limit. However, one can already pre- 
dict that if the perturbation does not couple any of the 
two eigenstates corresponding to K with the two eigen- 
states corresponding to K , the system will remain gap- 
less at all orders for e = 0. Indeed, in this case, the 



single-cone approximation proposed in Refs. [6—8] can be 
made safely, leading to a finite gap only when a scalar 
as well as a vector potential are present. Thus, to open 
the gap, one must have a perturbing potential that cou- 
ples these two twofold-degenerate subspaces. Since this 
potential has, by construction, the same periodicity as 
C(p, q), this condition requires the existence of a recipro- 
cal lattice vector associated to 61 and b 2 , which equals 
K' — K. It is then straightforward to show that this 
condition is strictly equivalent to 



1 



= mod 3, 



(6) 



where v is the vortex density defined in Eq. (4). Dirac 
states then have a momentum k — mod (6 1 , b 2 ) where 
b* and b* 2 are primitive vectors of the reciprocal lattice 
associated to foi and b 2 (b*.bj = 2TtSi t j). Let us under- 
line that this is a necessary condition that might not be 
sufficient to open a gap but, as we shall see, it is. 

A naive first-order degenerate perturbation theory 
consists in considering the subspace spanned by the 
four Dirac states. There, one gets a finite gap 
A(fe = 0) = 2-k<\>^/v. However, it is clear that condition 
(6) together with the similar conic dispersions near K 
and K implies that states in the vicinity of the Dirac 
cones are also coupled by the perturbation and one must 
look for fe 7^ states that may have a lower gap. Of 
course, the corresponding subspace depends directly on 
the vector potential. For the gauge choice described af- 
ter Eq. (4), one finds that the state that has the low- 
est positive energy is found for feo = f b* 2 . We checked 
by exact diagonalizations that this remarkable result is 
valid for any flux <fi for the configurations C(p,q) with 
v < 1/12. The state with the lowest positive energy is 
therefore expected to always be found in this sector [11]. 
Of course, the corresponding energy may be degenerate 
and may also be found for other momenta, as is the case 
for = 1/2. 



1/v p q A/(7T0) 2 

3 11 1/3 

9 3 5/21 

12 2 2 1/6 

21 4 1 0.077586 

27 3 3 0.061324 

36 6 11/130 



TABLE I. Gap A, at order (j> 2 , for the first values of lju 
satisfying Eq. (6). For 1/v — 21 and 27, the gap cannot be 
expressed as a simple fraction and we only give the first digits 
obtained numerically. 

At first order in <fi, one gets A(feo) = so that the 
low-energy effect of the perturbation is simply to shift 
the Dirac cones [12] (without renormalizing the Fermi 
velocity at s = 0). To go beyond, one has to consider 



the second-order degenerate perturbation theory in the 
k = kg subspace. Such an analysis involves the compu- 
tation of matrix elements of the perturbation between 
all states belonging to this sector which, for arbitrary p 
and q, is not an easy task. The expression of the gap for 
the first fillings satisfying condition (6) is given in Table 
I. Although in general it is difficult to get a simple ex- 
pression of A, one can derive exact formulas for q = 
(p being a multiple of 3) that allow one to (numerically) 
investigate large unit cell systems that would be out of 



reach with exact diagonalizations. From now on, we will 
mainly focus on this subset of configurations for which 
the gap reads 



with 



A(p,g = 0) 

W) 2 



= C P ~\Bl 



•S\p 



(7) 



p-l 2p-l 



A>= iEE l^){ 3 + 4cos K? + l)]+ 2c °s[2-(7 " §)]}{l + «»( 2 f)+«»[£(2»-»0]}. (8) 



B. 



2« 4 ^-^ ^-^ s 2 (m,n) 

1 n=0 m=l y ' ' 



p-l 2p-l 



p = ^i 2^ 2^ 



£(p, m) 



2r,4 z_^ z_^ s^(m.n) 

1 n=0 m=l v ' ' 
p-12p-l , . 

p „4 Z^ Z^ e 2 (m,n) 

y n=0 m=l v ' ; 



{ 3 + 4cos [tt(* + §)]+ 2o os[27r(* -|)]}{ sin (^) + sin [f(2n-m)]}, (9) 

cos (^a) - cos (S)] {i + cos (?fi) + cos [f (2n - m)] }, (10) 



where the sum over to is performed over odd in- 
tegers only. For convenience, we also introduced 



1.-.0 



e (m,n) 



£(p, to) 



:(k 



ay* 



|a|) [see Eq. (5)], and 



/ 37rm \ 
V 2p / 



2p 

l-(-l) p/3 



if to 7^ mod 
otherwise. 



(11) 



In the large-p limit, it is clear that A vanishes since one 
has to recover the spectral properties of the zero-flux 
problem. To analyze this infinitely diluted vortex limit, 
we computed the gap using Eq. (7) up to p = 20000. A 
close inspection of A p , B p , and C p led us to conjecture 
that the gap vanishes as A/0 2 ~ vYn.v~ l in the large- 
p = l/y/v limit. A convincing check of this result is, how- 
ever, displayed in Fig. 2. A natural question that arises at 
this stage concerns the behavior of the gap away from the 
perturbative regime analyzed up to now. To investigate 
arbitrary fluxes, one must diagonalize H numerically but 
the main advantage is that one only has to consider the 
subspace corresponding to k = k [] L] where the lowest- 
positive energy state lies. However, for arbitrary fluxes, 
one is restricted to small values of p since the number of 
sites per unit cell is Ap 2 and we need the full spectrum 
of the k = k subspace. In Fig. 3, we display the be- 
havior of the gap as a function of <fi and for q = and 
p = 3, 6, . . . , 51. As can be seen, A is a monotonously 
decreasing (increasing) function of p (of <fi in the interval 
[0,1/2]). We have also observed that the way A vani- 
shes when p increases depends on <f>. However, the lack 
of large-p data prevents to perform a sound analysis of 
these behaviors. 




FIG. 2. (Color online) Behavior of p 2 A/(f> 2 as a function of 
hip (for q — 0) in the small-0 limit. Exact results are obtained 
from Eq. (7) and the full line is a linear fit in good agreement 
with the conjecture discussed in the text. 



As already observed in the small-</> limit (see Table I), 
the maximum value of the gap is obtained for the largest 
vortex density satisfying Eq. (6), i.e., v = 1/3, but it is 
also obtained for the largest possible flux, i.e., cj> = 1/2. 
Denoting x* , the smallest positive root of the following 
polynomial 

P(x) = x 6 - 18 x 5 + 117 x A - 340 x 3 + 428 x 2 - 176 x + 16, 

(12) 
one gets 

A(p=l,q = l,<f>= 1/2) = 2s/x* ~ 0.70884. (13) 

Surprisingly enough, the problem considered here for 




51 



FIG. 3. (Color online) A as a function of 4> an d p (for q — 0). 
The maximum is reached for p = 3 and <j> — 1/2 where 
A ~ 0.611132. 



<p = 1/2 is directly connected to the spin 1/2 model intro- 
duced by Kitaev in 2006 [9], Using a Majorana fermionic 
representation of Pauli matrices, the Kitaev model can 
be mapped onto a free-fermion model in a 1^ gauge field. 
As a consequence, the value of the (effective) flux in each 
elementary plaquette is restricted to <p = or 1/2. This 
correspondence allowed Kitaev to identify the vortex con- 
figuration where the ground state of his system (Fermi 
sea at half-filling in the present electron language) lies. 
Indeed, as suggested early on in the flux-phase frame- 
work [13-15], the lowest energy at half- filling is obtained 
for = 0. In this problem, the more general question 
was as follows: for a given electron density, what is the 
flux density (and the flux pattern) that minimizes the en- 
ergy ? Although the answer has been provided by Lieb 
[16] for the special case of half-filling, exact results are 
still missing for arbitrary electron density. 

The present study raises a complementary question: 
given a flux density, what is the flux pattern that max- 
imizes the gap ? Undoubtedly, this question is even 
more difficult and the answer likely depends on the elec- 
tron density. For the Kitaev model (half-filling and 
<f> = 0, 1/2), we investigated several periodic configura- 
tions corresponding to fixed flux density v satisfying (6), 
and we are led to conjecture that the flux pattern maxi- 
mizing the gap is always C(p, q). Note that this flux pat- 
tern also minimizes the energy. One way to understand 
this result is to argue that the vortex-vortex interaction 
for (f> = 1/2 is repulsive so that it seems natural to find 
a triangular (Abrikosov-like) superlattice as an optimal 
pattern. However, it would be valuable to prove this 
result rigorously as well as finding gapped flux configu- 
rations for arbitrary v. It would also be worth adding 
further hopping processes as discussed in Ref. [17] that 
may give rise to a nontrivial insulator. Such consider- 
ations are clearly beyond the scope of the present work 
but we hope to have underlined that interesting phenom- 
ena may occur for nontrivial vortex configurations in the 
honeycomb lattice (see also Ref. [18] for related studies 
of the Kitaev model) . An obvious consequence of our re- 
sults for the Kitaev model is that there must be a finite 



gapped region around the point where a gap is induced 
by the vortex superlattice C(p, q). This is due to the fact 
that an insulator, as the one considered here, is robust 
to small deformations (for example, anisotropics in the 
hopping elements). 

One must also wonder how to observe this metal- 
insulator transition induced by a vortex superlattice in 
the honeycomb lattice. Obviously, the main difficulty is 
the realization of the superlattice with the ad hoc pa- 
rameters. The most realistic choice would be a flux 
= 1/2 for which vortices and antivortices are equiv- 
alent, so that one can use the exact correspondence to 
the Kitaev model for which many experimental propos- 
als exist [19-22]. Otherwise, in the context of graphene, a 
type-II superconductor might be used in the mixed state 
where the Abrikosov vortex lattice is found. Then, given 
6 = 1/2, one could think about gluing a graphene sheet 
on top of the superconductor. However, one faces the 
problem that the vortex core is much larger than a sin- 
gle elementary plaquette. Therefore, one has to investi- 
gate the gap-opening problem in the presence of extended 
though localized flux spots. Another appealing approach 
would be to consider optical flux lattices recently sug- 
gested in Refs. [23, 24] that seem especially adapted to 
our problem. 

Finally, given the occurrence of Dirac points in 
many experimental devices (see Ref. [ ] for a recent 
discussion) , we hope that the present work will motivate 
further investigations concerning the fate of these 
singularities in the presence of a vortex superlattice. 
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